home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
FishMarket 1.0
/
FishMarket v1.0.iso
/
fishies
/
301-325
/
disk_322
/
gwin
/
examples
/
geomap.c
< prev
next >
Wrap
C/C++ Source or Header
|
1992-05-06
|
13KB
|
594 lines
#include "gwin.user.h"
#include <math.h>
double myatan2();
float t0[3],t1[3],t2[3],ti0[3],ti1[3],ti2[3];
double pi,pi1, pi2, pi3;
float xmerc();
float lats[500],longs[500];
int ortho;
int merc;
int setmercator(),setortho(),plot_grid(),interrupt_display();
int new_center(),redisplay(),setortho2();
float x,y;
int event;
char key;
int interrupt = 0;
int display_underway = 0;
float lat_cen,long_cen;
int ortho2;
main()
{
float xlat,xlong;
char line[30];
pi1 = 3.1415926535897943/180.;
pi2 = 180./3.1415926535897943;
pi3 = 3.1415926535897943/2.0;
ortho = 1;
ortho2 = 0;
merc = 0;
printf("\n\nEnter lat_cen, long_cen: ");
scanf("%f %f",&lat_cen,&long_cen);
ustart("high2",0.,640.,0.,400.);
if(merc ) uwindo(-180.,180.,-180.,180.);
if(ortho || ortho2) uwindo(-250.,250.,-160.,160.);
setup_menus();
make_rotation_matrices(lat_cen,long_cen);
make_geo_grid();
plot_grid();
while(1 == 1) {
ugrina(&x,&y,&event,&key);
latlong(x,y,&xlat,&xlong);
sprintf(line,"%6.2f %6.2f",xlat,xlong);
if (ortho || ortho2){
uprint(100.,150.,line);
}else{
uprint(100.,170.,line);
}
}
}
latlong(x1,y1,xlat,xlong)
float x1,y1,*xlat,*xlong;
{
float xxx,yyy,zzz;
double xtemp,ytemp,ztemp;
float temp1,temp2;
float ph1,th1;
char line[30];
if(ortho || ortho2){
yyy = x1/150.;
zzz = y1/150.;
temp1 = yyy*yyy;
temp2 = zzz*zzz;
if(temp1 + temp2 > 1.0)return(1);
xxx = sqrt((double)(1 - temp1 - temp2));
xtemp = (double)(ti0[0] * xxx + ti0[1] * yyy + ti0[2] * zzz);
ytemp = (double)(ti1[0] * xxx + ti1[1] * yyy + ti1[2] * zzz);
ztemp = (double)(ti2[0] * xxx + ti2[1] * yyy + ti2[2] * zzz);
if(fabs(xtemp) < 1e-15) xtemp = 1e-15;
*xlong = pi2*myatan2(ytemp,xtemp);
if(ztemp > 1.0) ztemp = 1.0;
if(ztemp < -1.0) ztemp = -1.0;
*xlat = pi2*(pi3 - acos(ztemp));
}
if(merc){
ph1 = pi1 * y1;
ph1 = 2.0 * atan(exp(ph1)) - pi3;
ph1 = pi3 - ph1;
th1 = pi1 * x1;
xxx = sin((double)ph1) * cos((double)th1);
yyy = sin((double)ph1) * sin((double)th1);
zzz = cos((double)ph1);
xtemp = (double)(ti0[0] * xxx + ti0[1] * yyy + ti0[2] * zzz);
ytemp = (double)(ti1[0] * xxx + ti1[1] * yyy + ti1[2] * zzz);
ztemp = (double)(ti2[0] * xxx + ti2[1] * yyy + ti2[2] * zzz);
if(fabs(xtemp) < 1e-15) xtemp = 1e-15;
*xlong = pi2*myatan2(ytemp,xtemp);
if(ztemp > 1.0) ztemp = 1.0;
if(ztemp < -1.0) ztemp = -1.0;
*xlat = pi2*(pi3 - acos(ztemp));
}
}
make_rotation_matrices(lat_cen,long_cen)
float lat_cen,long_cen;
{
/* Rotation matrices to rotate lat/long to new lat/long */
/* so that (lat_cen,long_cen) becomes (0,0) */
float sin_theta, cos_theta, theta_rot, sin_phi, cos_phi, phi_rot;
theta_rot = long_cen * pi1;
phi_rot = lat_cen * pi1;
sin_theta = sin(theta_rot);
cos_theta = cos(theta_rot);
sin_phi = sin(phi_rot);
cos_phi = cos(phi_rot);
t0[0] = cos_theta * cos_phi;
t0[1] = sin_theta * cos_phi;
t0[2] = sin_phi;
t1[0] = -sin_theta;
t1[1] = cos_theta;
t1[2] = 0.0;
t2[0] = -cos_theta * sin_phi;
t2[1] = -sin_theta * sin_phi;
t2[2] = cos_phi;
ti0[0] = cos_theta * cos_phi;
ti0[1] = -sin_theta;
ti0[2] = -cos_theta * sin_phi;
ti1[0] = sin_theta * cos_phi;
ti1[1] = cos_theta;
ti1[2] = -sin_theta * sin_phi;
ti2[0] = sin_phi;
ti2[1] = 0.0;
ti2[2] = cos_phi;
}
mercator(long_in,lat_in,x_out,y_out)
float lat_in,long_in,*y_out,*x_out;
{
float xla,xlo,x[3],xr[3],sth,cth,sph,cph;
double th1,ph1;
xla = lat_in;
xlo = long_in;
th1 = pi1 * xlo;
ph1 = pi3 - xla * pi1;
sth = sin(th1);
cth = cos(th1);
sph = sin(ph1);
cph = cos(ph1);
x[0] = sph * cth;
x[1] = sph * sth;
x[2] = cph;
xr[0] = t0[0] * x[0] + t0[1] * x[1] + t0[2] * x[2];
xr[1] = t1[0] * x[0] + t1[1] * x[1] + t1[2] * x[2];
xr[2] = t2[0] * x[0] + t2[1] * x[1] + t2[2] * x[2];
*y_out = xmerc( pi2 * asin(xr[2]));
*x_out = pi2 * myatan2(xr[1],xr[0]);
}
float xmerc(xlat)
float(xlat);
{
double a1;
if (fabs(xlat) < 89.9999){
if(xlat >= 0.0){
a1 = pi1 * (45.0 + .5 * xlat);
return((float) pi2 * log(tan(a1)));
} else {
a1 = -pi1 * (-45.0 + .5 * xlat);
return((float) -pi2 * log(tan(a1)));
}
} else {
if(xlat < 0) return(-750.);
return(750.);
}
}
orthographic(long_in,lat_in,x_out,y_out,frontside)
float lat_in,long_in,*y_out,*x_out;
int *frontside;
{
float xla,xlo,x[3],xr[3],sth,cth,sph,cph;
double th1,ph1;
xla = lat_in;
xlo = long_in;
th1 = pi1 * xlo;
ph1 = pi3 - xla * pi1;
sth = sin(th1);
cth = cos(th1);
sph = sin(ph1);
cph = cos(ph1);
x[0] = sph * cth;
x[1] = sph * sth;
x[2] = cph;
xr[0] = t0[0] * x[0] + t0[1] * x[1] + t0[2] * x[2];
xr[1] = t1[0] * x[0] + t1[1] * x[1] + t1[2] * x[2];
xr[2] = t2[0] * x[0] + t2[1] * x[1] + t2[2] * x[2];
*y_out = 150.0 * xr[2];
*x_out = 150.0 * xr[1];
*frontside = 0;
if(xr[0]>0.0)*frontside = 1;
}
make_geo_grid()
{
int i;
float x;
i = 0;
for(x = -90.0;x < 91.0 ;x += 10.0){
lats[i++] = x;
}
lats[0] = -89.999;
lats[i-1] = 89.999;
i = 0;
for(x = -180.0; x < 181.0; x += 10.0){
longs[i++] = x;
}
longs[0] = -179.999;
longs[i-1] = 179.999;
}
plot_grid()
{
if(merc){
mercplot();
return(0);
}
if(ortho){
orthoplot();
return(0);
}
if(ortho2){
ortho2plot();
return(0);
}
}
mercplot()
{
int i,j;
float x,y,xold;
int event,key;
display_underway = 1;
interrupt = 0;
uerase();
upset("colo",1.0);
for(i=0;i<19;i++){
mercator(longs[0],lats[i],&x,&y);
umove(x,y);
xold = x;
for(j=1;j<37;j++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
mercator(longs[j],lats[i],&x,&y);
if(fabs(xold-x) < 100.0)udraw(x,y);
umove(x,y);
xold = x;
}
}
for(j=0;j<37;j++){
mercator(longs[j],lats[0],&x,&y);
umove(x,y);
xold = x;
for(i=1;i<19;i++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
mercator(longs[j],lats[i],&x,&y);
if(fabs(xold-x) < 100.0)udraw(x,y);
umove(x,y);
xold = x;
}
}
interrupted:
interrupt = 0;
display_underway = 0;
}
orthoplot()
{
int i,j;
float x,y,xold,yold;
int event,key;
int frontside;
display_underway = 1;
interrupt = 0;
uerase();
upset("colo",1.0);
for(i=0;i<19;i++){
orthographic(longs[0],lats[i],&x,&y,&frontside);
umove(x,y);
xold = x;
yold = y;
for(j=1;j<37;j++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
orthographic(longs[j],lats[i],&x,&y,&frontside);
if(fabs(xold-x) < 100.0 && frontside){
umove(xold,yold);
udraw(x,y);
}else{
umove(x,y);
}
xold = x;
yold = y;
}
}
for(j=0;j<37;j++){
orthographic(longs[j],lats[0],&x,&y,&frontside);
umove(x,y);
xold = x;
yold = y;
for(i=1;i<19;i++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
orthographic(longs[j],lats[i],&x,&y,&frontside);
if(fabs(xold-x) < 100.0 && frontside){
umove(xold,yold);
udraw(x,y);
}else{
umove(x,y);
}
xold = x;
yold = y;
}
}
interrupted:
interrupt = 0;
display_underway = 0;
}
ortho2plot()
{
int i,j;
float x,y,xold,yold;
int event,key;
int frontside;
frontside = 1;
display_underway = 1;
interrupt = 0;
uerase();
/* plot backsides first */
upset("colo",10.0);
for(i=0;i<19;i++){
orthographic(longs[0],lats[i],&x,&y,&frontside);
umove(x,y);
xold = x;
yold = y;
for(j=1;j<37;j++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
orthographic(longs[j],lats[i],&x,&y,&frontside);
if(fabs(xold-x) < 100.0 && !frontside){
umove(xold,yold);
udraw(x,y);
}else{
umove(x,y);
}
xold = x;
yold = y;
}
}
for(j=0;j<37;j++){
orthographic(longs[j],lats[0],&x,&y,&frontside);
umove(x,y);
xold = x;
yold = y;
for(i=1;i<19;i++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
orthographic(longs[j],lats[i],&x,&y,&frontside);
if(fabs(xold-x) < 100.0 && !frontside){
umove(xold,yold);
udraw(x,y);
}else{
umove(x,y);
}
xold = x;
yold = y;
}
}
/* now plot front sides */
upset("colo",1.0);
for(i=0;i<19;i++){
orthographic(longs[0],lats[i],&x,&y,&frontside);
umove(x,y);
xold = x;
yold = y;
for(j=1;j<37;j++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
orthographic(longs[j],lats[i],&x,&y,&frontside);
if(fabs(xold-x) < 100.0 && frontside){
umove(xold,yold);
udraw(x,y);
}else{
umove(x,y);
}
xold = x;
yold = y;
}
}
for(j=0;j<37;j++){
orthographic(longs[j],lats[0],&x,&y,&frontside);
umove(x,y);
xold = x;
yold = y;
for(i=1;i<19;i++){
ugrina(&x,&y,&event,&key);
if(interrupt) goto interrupted;
orthographic(longs[j],lats[i],&x,&y,&frontside);
if(fabs(xold-x) < 100.0 && frontside){
umove(xold,yold);
udraw(x,y);
}else{
umove(x,y);
}
xold = x;
yold = y;
}
}
interrupted:
interrupt = 0;
display_underway = 0;
}
setup_menus()
{
uamenu(1,0,0,"PROJECTION ",' ',0,(USHORT)MIDRAWN|MENUENABLED,7);
uamenu(1,1,0," ORTHOGRAPHIC ",'O',6,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED|CHECKIT,setortho);
uamenu(1,2,0," MERCATOR ",'M',5,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED,setmercator);
uamenu(1,3,0," ORTHO2 ",'2',3,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED,setortho2);
uamenu(2,0,0,"DISPLAY ",' ',0,(USHORT)MIDRAWN|MENUENABLED,7);
uamenu(2,1,0,"REDISPLAY ",'R',0,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED,redisplay);
uamenu(2,2,0,"INTERRUPT ",'I',0,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED,interrupt_display);
uamenu(2,3,0,"NEW CENTER ",'C',0,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED,new_center);
}
redisplay()
{
plot_grid();
interrupt=1; /* cancel previous invocation if any */
}
setortho()
{
uwindo(-250.,250.,-160.,160.);
ortho = 1;
merc = 0;
ortho2 = 0;
uamenu(1,1,0," ORTHOGRAPHIC ",'O',6,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED|CHECKIT|CHECKED,setortho);
redisplay();
}
setmercator()
{
uwindo(-180.,180.,-180.,180.);
merc = 1;
ortho = 0;
ortho2 = 0;
uamenu(1,2,0," MERCATOR ",'M',5,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED|CHECKIT|CHECKED,setmercator);
redisplay();
}
setortho2()
{
uwindo(-250.,250.,-160.,160.);
ortho2 = 1;
ortho = 0;
merc = 0;
uamenu(1,3,0," ORTHO2 ",'2',3,(USHORT)MIDRAWN|ITEMTEXT|HIGHCOMP
|COMMSEQ|ITEMENABLED|CHECKIT|CHECKED,setortho2);
redisplay();
}
interrupt_display()
{
interrupt = 1;
}
new_center()
{
char coords[255];
uerase();
uwindo(0.,100.,0.,100.);
ugetstring(10.,75.,70.,"Enter latitude and longitude: ",coords);
while (sscanf(coords,"%f %f",&lat_cen,&long_cen) < 2){
ugetstring(10.,75.,70.,
"Error, try again. Enter latitude and longitude: ",coords);
}
make_rotation_matrices(lat_cen,long_cen);
if(merc)setmercator();
if(ortho)setortho();
if(ortho2)setortho2();
}
double myatan2(y,x)
double x,y;
{
double temp;
double pi = 3.1415926535897943;
double piov2 = 0.5 * 3.1415926535897943;
if(fabs(x) > 0.0){
temp = atan((double)y/x);
}else{
if(y < 0.0){
temp = -piov2;
}else{
temp = piov2;
if(fabs(x) < 1e-15 && fabs(y) <1e-15) temp = 0.0;
}
}
if(x < 0.0 ){
if(y >= 0.0) {
temp = pi + temp;
}else{
temp = -pi + temp;
}
}
return(temp);
}